#This is to get:
#Figure 1

rm(list = ls())
library("rgdal")
library("ggplot2")
library("dplyr")
library("raster")
library("RColorBrewer")
library("lemon")
library("sp")
library("ggsn")
library("RColorBrewer")
library("lemon")

setwd("/Users/bgpopescu/Dropbox/Legacies_Central_Europe/")
setwd("C:/Users/bogdanp/Dropbox/Legacies_Central_Europe/")

#Step1: Reading CVS file with no. years under the Ottoman Empire
ottoman_years<-read.csv(file="./data/ottoman_years_nuts.csv")

#Step2: Reading shapefile for NUTS polygons (counties)
nuts <- readOGR(dsn="./data/shape_files.gdb",
                layer="nuts_NONEU_and_EU", use_iconv = TRUE, encoding = "UTF-8")

#Step3: Reading shapefile for countries
all_coun <- readOGR(dsn="./data/shape_files.gdb",
                    layer="all_country_boundaries", use_iconv = TRUE, encoding = "UTF-8")

#Step4: Preparing files for mapping
nuts_fort <- fortify(nuts) 
nuts@data$id <- rownames(nuts@data)
mapa.df <- right_join(nuts_fort,nuts@data, by="id")
mapa.df$MERGE_ID<-as.numeric(as.character(mapa.df$MERGE_ID))
all_countries<- fortify(all_coun)

#Step5: Merging the csv file (no. years under the Ottomans with the nuts polygons)
mapa.df2 <- left_join(mapa.df, ottoman_years, by=c("MERGE_ID"="MERGE_ID"))
mapa.df2$years_ottoman[is.na(mapa.df2$years_ottoman)]<-0
rm(mapa.df, ottoman_years)
rm(all_coun, nuts)

#Step6: Making the map
#Step6a: Defining no. classes
no_classes <- 5
labels <- c()

#Step6b: Manually defining categories
obj <- c(1, 100, 200, 300, 400, 500, 593)

#Step6c: Defining custom labels
labels <- c()
for(idx in 1:length(obj)){
  labels <- c(labels, paste0(round(obj[idx], 4), 
                             "-", 
                             round(obj[idx + 1], 4)))
}

#Step6d: Removing the last label to avoid something like "66.62 - NA"
labels <- labels[1:length(labels)-1]

#Step6e: Creating cutoff points
mapa.df2$brks <- cut(mapa.df2$years_ottoman, 
                    breaks = obj, 
                    labels = labels, 
                    include.lowest = T)


#Step6f: Defining colors
my.cols <- brewer.pal(6, "Blues")

#Step6g: Making mapp
map<-ggplot() + 
  geom_polygon(data = mapa.df2, colour = "grey60", aes(x=long, y=lat, group = group, fill=brks))+
    geom_path(data = all_countries, aes(x=long, y=lat, group = group), colour = "black", size = 0.5)+
    #Define zoom
    coord_cartesian(ylim=c(35.00, 57.284524), xlim=c(9.479832,46))+
    theme_bw()+
    labs(x = "Longitude (X)", y="Latitude (Y)")+
    ggtitle("Years under the Ottoman Empire")+
    theme(plot.title = element_text(hjust = 0.5))+
    scale_fill_manual("Interval", values = my.cols, guide = "legend") +
    theme()+
    theme(legend.position = c(1, 0),
          #Legend.position values should be between 0 and 1. c(0,0) corresponds to the "bottom left"
          #and c(1,1) corresponds to the "top right" position.
          legend.box.background = element_rect(fill='white'),
          legend.background = element_blank(),
          legend.text=element_text(size=12))+
  scalebar(x.min = 38, x.max = 43,
           #Above you mention how long the scalebar should be
           y.min = 57, y.max = 57.1,
           #Above you mention how thick the scalebar should be
           dist = 200, dist_unit = "km",
           transform = T, model = "WGS84",
           location = "bottomright",
           st.size = 4,
           #Above you have the font size of the numbers below the scalebar
           st.dist =3,
           #Above you have the distance between the bar and the text, as a proportion of the y axis.
           height=3)

map2<-reposition_legend(map, 'bottom left')
ggsave(map2,file="./Paper/graphs/figure1.jpg", height=19.42, width=20.2, units = "cm", dpi=300)




